Effective Sampling in the Configurational Space 
by the Multicanonical-Multioverlap Algorithm 
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We propose a new generalized-ensemble algorithm, which we refer to as the multicanonical- 

Omultioverlap algorithm. By utilizing a non-Boltzmann weight factor, this method realizes a random 
walk in the multi-dimensional, energy-overlap space and explores widely in the configurational space 
' including specific configurations, where the overlap of a configuration with respect to a reference 

state is a measure for structural similarity. We apply the multicanonical-multioverlap molecular 
dynamics method to a penta peptide, Met-enkephalin, in vacuum as a test system. We also apply 
the multicanonical and multioverlap molecular dynamics methods to this system for the purpose 
of comparisons. We see that the multicanonical-multioverlap molecular dynamics method realizes 

S effective sampling in the configurational space including specific configurations more than the other 
two methods. From the results of the multicanonical-multioverlap molecular dynamics simulation, 
furthermore, we obtain a new local-minimum state of the Met-enkephalin system. 

a 

5a ' PACS numbers: 02.70.-c,05.10.-a,61.43.Fs 



I. INTRODUCTION 



In order to understand the protein folding mechanisms, it is essential that we investigate free-energy landscapes of 
Q ■ protein systems. In experiments, however, to obtain the free-energy landscapes of protein systems is very difficult. 
Therefore, computer simulations are now widely used for problems of the protein folding mechanisms. In computer 
simulations, however, it is still difficult to get accurate free-energy landscapes of protein systems. Accordingly, many 
I ' efforts are devoted to develop efficient simulation algorithms. 
, In complex systems such as proteins, we must realize effective samplings in the configurational space. In usual 
Q\ ' canonical-ensemble simulations Q, 0, 0, 0, IH El , however, it is difficult to achieve this. This is because the usual 
, canonical-ensemble simulations tend to get trapped in a few of many local-minimum states. To overcome these 
difficulties, the generalized-ensemble algorithms have been proposed (for a review, see, for instance, Ref. 0). 

The multicanonical algorithm 0, 0, 0, [nj is one of the most well-known methods among the generalized-ensemble 
algorithms. In the multicanonical ensemble, the probability distribution of the potential energy is expressed by 
the product of the density of states and a non-Boltzmann weight factor, which we refer to as the multicanonical 
weight factor, and we have a flat probability distribution of the potential energy. Therefore, multicanonical-ensemble 
simulations realize free random walks in the potential-energy space and have effective samplings in the configurational 
space. The random walk in the potential-energy space allows one to calculate various thermodynamic quantities 
I ■ as functions of temperature for a wide temperature range from the results of a single simulation run by the single- 
histogram reweighting techniques |!2lll3| . This method is suitable to sample widely the configurational space, but not 
to have samplings that focus on a specific configuration because of the very nature of the algorithm. Consequently, 
it is difficult to obtain accurate free-energy landscape around specific configurations in multicanonical simulations, 
while we can obtain the free-energy landscape over wide areas. 

The multioverlap algorithm was recently proposed for Monte Carlo (MC) method 0] and molecular dynamics 
(MD) method [TEl flfjj l in order to investigate the stability of specific configurations and transition states among 
specific configurations. In the multioverlap ensemble, the probability distribution is expressed by the product of the 
density of states and a non-Boltzmann weight factor, which we refer to as the multioverlap weight factor, and we 
have a flat probability distribution in the overlap space. The method aims at achieving effective samplings that 
focus on specific configurations. Accordingly, we can obtain an accurate free-energy landscape around the specific 
configurations and estimate correctly transition states among the specific configurations. In this method, however, 
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we do not have a random walk in the potential-energy space, and hence, thermodynamic quantities can be obtained 
at temperatures only near the simulation temperatures. 

In this article, we propose a simulation method, which we refer to as the multicanonical-multioverlap method, to 
sample widely the configurational space and effectively the vicinity of specific conformations of a protein. We apply 
the multicanonical-multioverlap MD method to Mct-cnkcphalin in vacuum and test the effectiveness of the method 
by comparing the results with those of the multicanonical MD method and the multioverlap MD method. 

In Section [H] we summarize the formulation of the multicanonical-multioverlap algorithms. We present the details 
of the three simulations that we performed and their results in Section liTTl Section Hvl is devoted to conclusions. 

II. MULTICANONICAL-MULTIOVERLAP ALGORITHMS 

In this section we describe the multicanonical-multioverlap algorithms. In Sec. Ill Al we define the dihedral-angle 
distance |l4Lll7| which is a complementary quantity of the overlap. Realizing effective samplings in the energy-overlap 
space, we introduce the multicanonical-multioverlap weight factor in Sec. Ill Bl In Sec. Ill Dl we present the equations 
of motion to perform multicanonical-multioverlap MD simulations. We describe multicanonical-multioverlap MC 
method with the Metropolis criterion p| in Sec. Ill Cl In Sec. Ill El we explain the reweighting techniques [I3.ll3j . We 
can obtain appropriate physical quantities at any temperature in the natural ensemble by the reweighting techniques. 

A. Definition of dihedral-angle distance 

We introduce a dihedral- angle distance, which is a complementary quantity of the overlap, as a reaction coordinate. 
The dihedral-angle distance d with respect to a reference configuration is defined by 

1 " 

d = —J2 d ^ v ^- 

mr z — ' 

i=l 

Here, n is the total number of dihedral angles, Vi is the dihedral angle i, and w° is the dihedral angle i of the reference 
configuration. The distance d a (vi,vf) between two dihedral angles is given by 

da(«i,«°) = min(K - u?|,27r - \ Vi - »?|) . (2) 

If d = 0, from Eq. QJ, all dihedral angles are coincident with those of the reference configuration, and the two 
structures are identical. The dihedral-angle distance is thus an indicator of how similar the conformation is to the 
reference conformation. 

Sampling the vicinity of the reference configuration is equivalent to sampling around d — 0. In order to realize a 
sampling that focuses on the reference configuration, therefore, we just have to sample the neighborhood of d = 0. 
We then would like to realize an effective sampling, which covers widely the configurational space and focuses on the 
reference configuration. In other words, we want to sample effectively and widely the configurational space including 
near d = 0. 



B. Effective sampling in the energy-overlap space 

In the case of canonical ensemble at a constant temperature To, the probability distribution P c of potential energy 
E is represented by the product of the density of states n(E) and the Boltzmann weight factor W C (E; To): 

P C (£;T ) = n(E)W c (E;T Q ) 

= n(E)e~P° E , (3) 

where [3$ is given by (3q = l/fceTo (fee is the Boltzmann constant). In the multicanonical ensemble at a constant 
temperature T , by employing the non-Boltzmann weight factor W muca (E), which we refer to as the multicanonical 
weight factor, a uniform probability distribution of potential energy is obtained: 

P m uca{E) = n(E)W muca {E) 

= constant , (4) 
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where E muca (E) is the multicanonical potential energy. The multicanonical weight factor, or the multicanonical 
potential energy, is not a priori known and has to be determined by sho rt p reliminary simulations. There exist many 
methods for the determination of the weight factor (see, e.g., Refs. [IH Ha, [S3, E3, [23, [2^, E3, H^, [2^|). Equation © 
implies that Therefore, multicanonical simulations realize a free random walk in the potential-energy space and are 
able to sample effectively the configurational space. In this method, however, it is difficult to sample the vicinity of 
specific configurations because of the very nature of the algorithm. 

The multioverlap algorithm, which is developed by generalizing the multicanonical algorithm, is suitable to have 
samplings that focus on specific configurations. In the multioverlap ensemble at a constant temperature To, the 
probability distribution of dihedral-angle distances is defined by 



Pmuov(di, ■ ■ ■ , djv) — / dE n(E; d±, ■ ■ ■ , djsr)W m uov(Ei d±, ■ ■ ■ , div) 

= constant , (5) 

where di is the dihedral- angle distance with respect to reference configuration i (i = 1, • • • , N), n(E; d\, ■ • ■ , djv) is 
the density of states, W muov (E] di, ■ ■ ■ , djv) is the multioverlap weight factor, and f(d±, • • • , djv) is the "dimensionlcss 
free energy" . The dimensionless free energy is not a priori known and has to be determined by short preliminary 
simulations. This method performs a random walk in the iV-dimensional dihedral-angle-distance space, in which the 
simulation visits the N reference configurations often. We can thus obtain accurate information about transition states 
among these N states. Because the multioverlap method does not realize a free random walk in the potential-energy 
space, it is difficult to sample widely the potential-energy space. 

We want simulations to have effective and wide samplings including near d = in the configurational space. 
Therefore, we consider to carry out a simulation that performs a random walk both in the potential-energy space and 
in the dihedral-angle-distance space (energy-overlap space). In other words, the simulation needs to have a constant 
probability distribution in the energy-overlap space. In analogy with the multicanonical ensemble in Eq. 10} or the 
multioverlap ensemble in Eq. J^J), by employing the non-Boltzmann weight factor W m co(E, di, ■ ■ ■ , d/v), which we refer 
to as the multicanonical-multioverlap weight factor, a uniform probability distribution with respect to the potential 
energy and dihedral-angle distances is obtained: 

P mco (E,di, ■ ■ ■ ,d N ) = n(E, di, ■ ■ ■ ,d N )W mco (E, d\, ■ ■ ■ ,d N ) 

= constant . (6) 

The multicanonical-multioverlap weight factor is not a priori known and again has to be determined by short prelim- 
inary simulations. The multicanonical-multioverlap weight factor W m co(E, di, ■ ■ ■ , d/v) at a constant temperature To 
can be written as 

W mco (E, d u ■ ■ ■ , d N ) = e-/WWJMi,- ,<M , (7) 

where E mco (E, di, ■■ ■ , djv) is the multicanonical-multioverlap potential energy. We remark that by definition the 
multicanonical weight factor in Eq. Q and the multicanonical-multioverlap weight factor in Eq. © are independent 
of temperature, whereas the multioverlap weight factor in Eq. J5J depends on temperature. 



C. Monte Carlo methods in the multicanonical-multioverlap ensemble 

Canonical MC simulations are performed with the Metropolis criterion 0] . In the Metropolis criterion, the transition 
probability from state x with potential energy E to state x' with potential energy E' is given by 

t >> - / 1 ' fOT AE - ' (R\ 

W(X X > ~ \ exp(-/3 A£) , for AE > , W 



where 



AE = E' - E . (9) 



In multicanonical-multioverlap MC simulations, the transition probability from state x with potential energy E 
and dihedral- angle distances di, • ■ • ,d^ to state x' with potential energy E' and dihedral-angle distances g? 1; ■ • ■ , d' N 
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is correspondingly given by 



w(x -+x'}-i 1 ' f ° r AE,nco - ' (10) 

1 J_i exp(-/3 A£ mco ) , for AE mco > , [W > 



where E mco is the multicanonical-multioverlap potential energy in Eq. J7J) and 

Ai? mco = E mco (E , dx, ■ • • , djv) — E mco (E, d±, ■ ■ ■ , <%) • (11) 

D. Equations of motion in the multicanonical-multioverlap ensemble 

Solving regular Newton's equations of motion leads to the microcanonical ensemble. There are several methods 
to realize the canonical ensemble by the MD simulation (see, for example, Refs. Q, 0, 0, 0, S El)- Here, we just 
consider one of these methods, namely the Gaussian constraint method |2j,y| (the realization by other methods is also 
straightforward). In the Gaussian constraint method, the following equations of motion with Gaussian thermostat 
are solved: 

= d qi = Pi 

dt rm ' (12) 
Pi =Fi- CcP l , 

where mi, q i , and p i are the mass, coordinate vector, and momentum vector of atom i. The force Fi acting on atom 
i is given by 

where E is the potential energy. The coefficient £ c is chosen so as to guarantee that the total kinetic energy is constant: 

£*w< 

Co = -r- • (14) 

i 

Correspondingly, the molecular dynamics algorithm in the multicanonical-multioverlap ensemble naturally follows 
from Eq. Q (see Refs. [HIHl for the case of multicanonical MD). The multicanonical-multioverlap MD simulation 
is carried out by solving the following modified equations of motion with Gaussian thermostat: 

■ _dqj _ Pj 

dt rrii (15) 

fi % \>mcol>i ■ 

The 'force' F™ co acting on atom i is calculated from 



7i//; co 



dE n 



where E mco is the multicanonical-multioverlap potential energy in Eq. J7J|. The coefficient £ mco is defined by 



(16) 



2m,i 
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E. Reweighting techniques 

The results of the multicanonical-multioverlap production run can be analyzed by the reweighting techniques. 
Suppose that we have determined the multicanonical-multioverlap potential energy E mco at a constant temperature 
T and that we have made a production run at this temperature. The expectation value of a physical quantity A at 
any temperature T is calculated from 



2J A(E,di,- ■ ■ ,d N )n(E,di,- ■ ■ ,d N )e 



-fiE 



< A> T = E ' dl '-' dN — , (18) 

n{E,d lt "- ,d N )e-? E 

E,d%,--- ,d]\ 



where the best estimate of the density of states is given by the single-histogram reweighting techniques 0, 0] (see 
Eq. ©): 

/ rp j 7 \ N mco (E, di, ■ ■ ■ , d N ) 

n{E,d lr -- ,d N ) = — — — — , (19) 

W mco {E,di, ■ ■ ■ ,d N ) 

and N mco (E, d\ 1 ■ ■ ■ , d/v) is the histogram of the probability distribution that was obtained by the multicanonical- 
multioverlap production run. By substituting Eqs. J7J) and <|19|) into Eq. I|18l) . we have 

]T A(E,d u --- ,d N )N mco (E,d lr -- ,d N y° E -° ^,d u ...,d N) -pE 

< A > T = E ^--- dN . (20) 

£ N m UE,d 1 ,---,d N )e^ E ^ E ' d ----^ E 

E,di,--- ,djy 

We can also calculate the free energy (or, the potential of mean force) with appropriate reaction coordinates. For 
example, the free energy F(^i,^2]T) with reaction coordinates £i and £2 at temperature T is given by 

F(^fr,T) = -fc B TlnP c (£i,6;T) , (21) 
where P c (£i, £,2',T) is the reweighted canonical probability distribution of £1 and £2 and given by (see Eq. J20J0 



N mco (£i,fr,E,di,- ■ ■ ,d N )e 



/3 E mco (E,d lt - 4n)-0e 



pic C .rp\ _ E,d u --- ,d N . . 

ilfy-E.dt,--- ,d N 

and N meo (£i, £2] E, d\, ■ •■ ,dw) is the histogram of the probability distribution that was obtained from the 
multicanonical-multioverlap production run. 



III. APPLICATION TO MET-ENKEPHALIN IN GAS PHASE 



A. Computational details 

In order to demonstrate the effectiveness of the multicanonical-multioverlap MD method, we compare a 
multicanonical-multioverlap MD simulation with multicanonical and multioverlap MD simulations. We apply the 
three simulations to the system of Met-enkephalin in vacuum. Met-enkephalin is one of the simplest peptides and has 
the amino-acid sequence Tyr-Gly-Gly-Phe-Met. This peptide is often adopted as a test system in biomolecular simu- 
lations. In our simulations the N-terminus and the C-terminus were blocked with the acetyl group and the N-methyl 
group, respectively. This is because we wanted the total charge of the Met-enkephalin system to be neutral. The 
force field that we adopted is the CHARMM param 22 parameter set j2]j . Leap-frog algorithm [2^| was employed for 
the numerical integration and the time step was taken to be 0.2 fs. The reason for using such a small time step is to 
perform simulations with high accuracy. 

In the multicanonical-multioverlap MD simulation, we must have a reference conformation. Therefore, we adopted 
the conformation in Fig. ^ as the reference conformation and set N = 1 in Eqs. © and Q. This conformation was 
obtained by the simulated annealing MD method [29j as follows: During the simulated annealing run, the temperature 
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was decreased linearly from 1000 K to 100 K with an increment of 50 K, and the canonical MD simulations were 
performed for 200 ps at each temperature. This simulated annealing MD run was repeated 10 times with different 
initial random numbers. The obtained final conformations were further minimized by the conjugate gradient method, 
and we finally got two conformations from the results of these minimizations. One of the two conformations was the 
reference conformation in Fig. ^ From our previous results of the multioverlap MD simulations |15l Il6| , we see that 
this reference conformation is one of the conformations in a local-minimum state and another conformation obtained 
by the results of these minimizations is the global-minimum state with the CHARMM param 22 parameter set. 

The backbone dihedral angles are of three types: the rotation angle around the N — C" bond of the backbone (</>), 
that around the C Q — C bond (ip), and that around the peptide bond C — N (w). Our multicanonical-multioverlap MD 
simulation was performed using the all-atom model, but we used only cf> and ip angles in the definition of the dihedral- 
angle distances in Eq. QJ. This is because the dihedral angles of the backbone u> have almost the fixed value of 180° 
for the peptide bond C — N. Furthermore, by using only the backbone dihedral angles (and not side-chain dihedral 
angles) as the elements of the dihedral-angle distances, we focused on the backbone structures of Met-enkephalin. In 
Eq. (JJJ, consequently, the number n of the elements of the dihedral-angle distances is 10 because Met-enkephalin has 
five pairs of <f> and ip. 

For the purpose of comparisons, we also performed a multicanonical MD simulation and a multioverlap MD sim- 
ulation for 9 ns at To = 300 K. We determined the multicanonical weight factor in Eq. Q as follows. We carried 
out canonical MD simulations at eight temperatures between 300 K and 1000 K with equal increment of 100 K and 
obtained ensemble-averages of the potential energy at each temperature. From the ensemble-averages of the potential 
energy, we calculated the derivative of the multicanonical potential energy 0j H3 : 



dE muca (E) 



T ° (23) 



T{E av 



dE 
with 

Eave = ( E ) T (E avc ) ■ ( 24 ) 

We integrated the derivative of the multicanonical potential energy and obtained the multicanonical weight factor. 
We adopted a random-coil conformation for the initial conformation of the multicanonical MD simulation production 
run. 

In the multioverlap MD simulation, we employed the 2-dimensional version of this method. In other words, we 
used two reference conformations in the multioverlap MD simulation and TV = 2 in Eq. (0. One of the two reference 
conformations was the conformation in Fig. ^ and the other one is shown in Fig. [21 This reference conformation in 
Fig. [His one of the conformations in the global-minimum state with the CHARMM param 22 par ameter set (see Ref. 
(la. ITql. We determined the multioverlap weight factor in Eq. JSJ by the following process [Tsl IT^|. Suppose that 
we have the dimensionless free energy / = /w in the Ith iteration of the short multioverlap MD simulation. In the 
/ + 1th iteration of the short multioverlap MD simulation, is calculated from 

f« +1 Xd 1 ,d 2 ) = f«\d u d 2 )-\ogN«\d u d 2 ) , (25) 

where d\ and d 2 is the dihedral-angle distance for Reference Conformation 1 (RC1) in Fig. ^ and Reference Confor- 
mation 2 (RC2) in Fig. [3 respectively. N^ L \di,d 2 ) in Eq. I|25|l is the histogram obtained from the results of the Ith 
iteration. For this calculation, the dihedral-angle distances were discretized with a bin size of 0.01. Moreover, we 
interpolated the dimensionless free energy by a polynomial, following the techniques that were introduced in Ref. |30| 
(see Eq. (94) there). The initial value was set as follows: 

fW(d 1 ,d 2 ) = 0. (26) 

In the first iteration, therefore, we performed a short usual canonical MD simulation. We stopped the iterations 
after seven short multioverlap MD simulations each for 1.8 ns and we obtained the multioverlap weight factor. The 
initial conformation for the multioverlap production run was a conformation equilibrated at 300 K by the canonical 
simulation. Since the multioverlap weight factor has temperature dependence, it is appropriate that we employ this 
initial conformation. 

The multicanonical-multioverlap MD simulation was carried out at To = 300 K. We first have to determine the 
multicanonical-multioverlap weight factor W mC o(E, d\) in Eq. JJJ to get a flat probability distribution in the energy- 
overlap space (E, di). For this purpose we used a similar procedure to that in Eq. (|25|) . Namely, suppose that we have 
Emco = Emco in the Zth iteration of the short multicanonical-multioverlap MD simulation. In the I + 1th iteration of 
the short multicanonical-multioverlap MD simulation, Emto is calculated from 

E%g>(E, di) = E&l (E, di) + k B T \ogM l \E, di) , (27) 
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where N^ l '(E,di) is the histogram obtained from the results of the Ith iteration. For this calculation, the potential 
energy and the dihedral-angle distance were discretized with a bin size of 1.0 kcal/mol and a bin size of 0.01, respec- 
tively. We also interpolated the multicanonical-multioverlap potential energy by the polynomial. The initial value 
was set as follows: 

E mC o(E,di) — E muca (E) , (28) 

where E muca (E) is the multicanonical weight factor that was determined as above. We then performed three iterations 
of the multicanonical-multioverlap MD simulations in Eq. (|15|) for 3 ns. The multicanonical-multioverlap weight 
factor E mco (E,di) was updated by Eq. Ij27(l after each multicanonical-multioverlap MD simulation. Finally, the 
multicanonical-multioverlap MD production run was then performed with this weight factor for 9 ns after equilibration 
of 1 ns. For the initial conformation of the multicanonical-multioverlap MD simulation production run, we also adopted 
a random-coil conformation. 



B. Comparisons of the three simulations 

We first compare the time series of the potential energy and the dihedral-angle distance obtained from the multi- 
canonical, multioverlap, and multicanonical-multioverlap MD simulations. Fig.[3]shows the time series of the potential 
energy of the three simulations. The multicanonical and multicanonical-multioverlap MD simulations cover widely the 
potential-energy space, as we can see in Figs.|3Ia) and^c). In other words, the two simulations realized free-random 
walks in the potential-energy space and sampled widely the conformational space. In the multioverlap MD simulation, 
however, we can sample only a narrow region in the potential-energy space as in Fig. UJb). Therefore, in contrast 
with the other two simulations, the multioverlap MD simulations are not suitable to sample widely the conformational 
space. 

In Fig.0]we show the time series of the dihedral-angle distance d\ with respect to Conformation RC1 in Fig.Q] When 
di = 0, the values of the backbone dihedral angles are completely coincident with those of RC1. In the multioverlap 
and multicanonical-multioverlap MD simulations, we see from Figs. Hfb) and 0f c) that the efficient samplings were 
realized in the neighborhood of d\ = 0. In other words, the multioverlap and multicanonical-multioverlap MD 
simulations could sample efficiently the vicinity of RC1. The multicanonical MD simulation sampled infrequently 
the neighborhood of d\ = 0, as we can see in Fig. 0Ja). Thus, it is difficult to sample specific conformations in 
multicanonical MD simulations. 

From Fig. |3{c) and Fig. 0Jc) we see that the multicanonical-multioverlap MD simulation sampled both widely in 
the conformational space and efficiently the vicinity of the reference conformation. Therefore, the multicanonical- 
multioverlap MD method has the advantages of both the multicanonical MD method and the multioverlap MD 
method. 

We now consider the probability distributions of configurations from the three simulations. In Fig. we show the 
raw data of the histograms with respect to the potential energy E and dihedral-angle-distance d%. From Fig. Ufa) 
it is apparent that the multicanonical MD simulation had only a partial sampling in the vicinity of RC1 (low d\ 
regions). Moreover, the multicanonical MD simulation did not sample widely in the conformational space at the 
low-energy region, although it had a wide sampling in the conformational space at the high-energy region. In the 
multioverlap MD simulation, on the other hand, the sampling that focuses on RC1 was realized (see Fig. EJb)). At 
the high-energy region, however, we could not sample at all in the multioverlap MD simulation. From Fig. [Sic) we see 
that the multicanonical-multioverlap MD simulation performed the effective sampling in the conformational space in 
comparison with the other two simulations. In fact, we could sample widely the conformational space at low-energy 
region as well as high-energy region and the vicinity of RC1 in the multicanonical-multioverlap MD simulation. 

C. Physical quantities calculated by the reweighting techniques 

We present various physical quantities calculated from the results of the three simulations by the reweighting 
techniques. The reweighting techniques for the multicanonical-multioverlap MD method were exp lained in Sec. Ill El 
The reweighting techniques for the other MD methods are accounted, for instance, in Refs. 

UEM In Fig. we show 

the average potential energy and specific heat calculated as functions of temperature by the reweighting techniques. 
The specific heat here is defined by 

^ = 1 d(E) T 
dT 

= 1 ((£ 2 ) T - (Ef T ) . (29) 
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From the Figure we see that the results from the multicanonical-multioverlap MD simulations well coincide with 
those from the multicanonical MD simulation. The results from the multioverlap MD simulation, however, are in 
agreement with those from the other two simulations only near T — 300 K. This is because the multioverlap method 
sampled only a narrow potential-energy range as shown in Fig. 0b) . Accordingly, we can obtain various physical 
quantities at any temperature in the multicanonical-multioverlap MD method as in the multicanonical MD method, 
although it is difficult to get such quantities at any temperature in the multioverlap MD method. 

Subsequently, we present the free-energy landscape with respect to various reaction coordinates. The free-energy 
landscape was calculated from Eq. (|21|) with appropriate reaction coordinates by the reweighting techniques. In Fig.0 
we show the free-energy landscape at T — 300 K obtained from the three simulations with respect to the potential 
energy E and dihedral- angle distance d\. Since the multicanonical MD simulation did not sample the vicinity of RC1 
and widely the conformational space at the low-energy region, the free-energy landscape at T = 300 K was obtained 
only in a narrow region away from RC1 as shown in Fig. 0a). In the multioverlap MD simulation, on the other 
hand, the free-energy landscape in Fig. 0b) including the neighborhood of RC1 was calculated. From Fig. 0c), the 
free-energy landscape obtained from the multicanonical-multioverlap MD simulation was described over a wide range 
and near RC1. Moreover, we could identify a local-minimum state, which is located around (E,di) = (75.0,0.32) in 
Fig. 0c), although we did not find it by the other two methods. This is because the multicanonical-multioverlap MD 
simulation samples widely the conformational space at low-energy region and effectively in the vicinity of RC1. 

In Fig. [S] we also show the free-energy landscape at T — 300 K obtained from the three simulations with respect to 
the potential energy and root-mean-square distance (RMSD) n. Here, the RMSD r x with respect to RC1 is defined 



where N is the number of atoms, {q®} are the coordinates of RC1, and the minimization is taken over the rigid 
translations and rigid rotations of the coordinates of the configuration {q^. In this article, we took into account 
only the backbone coordinates of Met-enkephalin in Eq. i|30[l . In Fig. [S] we can also see a local- minimum state, which 
is located around [E, r{) = (75.0,2.0), only from the results of the multicanonical-multioverlap MD simulation. By 
analyzing the conformations in the two local minima found in Figs. 0c) and|HJc), we believe that they correspond to 
the same state. This is because common conformations were found in both states. In Fig. [5]we show a representative 
conformation in this newly found local-minimum state. This conformation has a backbone hydrogen bond between 
hydrogen bond donor NH of Gly-2 and hydrogen bond acceptor CO of Met-5. The conformation also has a backbone 
hydrogen bond between CO of Try-1 and NH of Gly-3 but this hydrogen bond was frequently broken. 



In this article, we have proposed the multicanonical-multioverlap MD algorithm, which is useful to sample the 
conformational space widely and the vicinity of a reference conformation effectively. We applied this method to a 
penta-peptide system of Met-enkephalin in vacuum and compared the performance with the those of multicanonical 
and multioverlap MD methods. We showed the effectiveness of the multicanonical-multioverlap MD method over the 
multicanonical and multioverlap MD methods. The multicanonical MD simulation sampled widely the conformational 
space at high-energy region but not at low-energy region and did not have the sampling around the reference confor- 
mation. The multioverlap MD simulation could sample effectively the vicinity of the reference conformation. In the 
multioverlap MD simulation, however, we were not able to have the sampling in the high-energy region. On the other 
hand, the multicanonical-multioverlap MD simulation realized a free-random walk in the energy-overlap space and 
sampled the conformational space widely and the neighborhood of the reference conformation. Accordingly, we could 
obtain accurate free-energy landscape in the wide reaction-coordinate space including the vicinity of the reference 
conformation and discover a new local-minimum state. 

In the protein folding problem, the multicanonical-multioverlap method can be applied to deduce folding pathways 
in which the protein system has an intermediate state like a molten-globule state. This is because we can obtain free- 
energy landscapes, which include random-coil states, the native state, and the molten-globule state, from the results 
of multicanonical-multioverlap simulations. Furthermore, we can estimate transition states accurately between the 
native state (or denatured state) and the molten-globule state by employing the molten globule state as the reference 
conformation in multicanonical-multioverlap simulations. 



by 
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FIG. 1: Reference configuration that was used in the multicanonical-multioverlap MD simulation. The dotted lines denote the 
hydrogen bonds. The N-terminus and the C-terminus are on the right-hand side and on the left-hand side, respectively. The 
figure was created with RasMol [3l|. 




FIG. 2: The other reference configuration that was used in the multioverlap MD simulation. See also the caption of Fig. 
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FIG. 3: The time series of the potential energy E. (a) is the results from the multicanonical MD simulation, (b) is from the 
multioverlap MD simulation at Tb = 300 K, and (c) is from the multicanonical-multioverlap MD simulation. 
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FIG. 4: The time series of the dihedral-angle distance di. (a) is the results from the multicanonical MD simulation, (b) is 
from the multioverlap MD simulation at To = 300 K, and (c) is from the multicanonical-multioverlap MD simulation. 
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FIG. 5: The raw data of the probability distribution with respect to the potential energy E and dihedral-angle distance di. 
(a) is the results from the multicanonical MD simulation, (b) is from the multioverlap MD simulation at To = 300 K, and (c) 
is from the multicanonical-multioverlap MD simulation. 
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FIG. 6: (a) Average potential energy as a function of temperature, and (b) specific heat as a function of temperature. 
These results were calculated from the multicanonical MD simulation (blue square), the multioverlap MD simulation (red 
triangle), and the multicanonical-multioverlap MD simulation (black circle) by the reweighting techniques. The results from 
the multicanonical and multicanonical-multioverlap simulations are essentially identical. 
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FIG. 7: The free-energy landscapes with respect to the potential energy E and dihedral-angle distance di that were obtained 
from (a) the multicanonical MD simulation, (b) the multioverlap MD simulation at To = 300 K, and (c) the multicanonical- 
multioverlap MD simulation. 
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FIG. 8: The free-energy landscapes with respect to the potential energy E and RMSD ri that were obtained from (a) the 
multicanonical MD simulation, (b) the multioverlap MD simulation at To — 300 K, and (c) the multicanonical-multioverlap 
MD simulation. 




FIG. 9: (a) A typical structure in the new local-minimum state found in Figs. |7Jc) and|8"Jc). (a) and (b) correspond to the 
same conformation viewed from different angles. See also the caption of Fig. Q 



